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The cellular basis for the adaptive immune response during antigen recognition relies on a spe¬ 
cialized protein interface known as the immunological synapse (IS). Understanding the biophysical 
basis for protein patterning by deciphering the quantitative rules for their formation and motion 
is an important aspect of characterizing immune cell recognition and thence the rules for immune 
system activation. We propose a minimal mathematical model for the physical basis of membrane 
protein patterning in the IS, which encompass membrane mechanics, protein binding kinetics and 
motion, and fluid flow in the synaptic cleft. Our theory leads to simple predictions for the spatial 
and temporal scales of protein cluster formation, growth and arrest as a function of membrane stiff¬ 
ness, rigidity and kinetics of the adhesive proteins, and the fluid in the synaptic cleft. Numerical 
simulations complement these scaling laws by quantifying the nucleation, growth and stabilization 
of proteins domains on the size of the cell. Direct comparison with experiment shows that passive 
elastohydrodynamics and kinetics of protein binding in the synaptic cleft can describe the short- 
time formation and organization of protein clusters, without evoking any active processes in the 
cytoskeleton. Despite the apparent complexity of the process, our analysis highlights the role of just 
two dimensionless parameters that characterize the spatial and temporal evolution of the protein 
pattern: a ratio of membrane elasticity to protein elasticity, and the ratio of a hydrodynamic time 
scale for fluid flow relative to the protein binding rate, and we present a simple phase diagram that 
encompasses the variety of patterns that can arise. 


INTRODUCTION 

Recognition of self or non-self is essential for an effective and functional adaptive immune response. The main 
players in this process are immune cells (T-lymphocyte cells (T-cells) B-cells, natural killer (NK) cells [1] and 

phagocytes 13111 that are constantly on the move scanning surfaces for antigenic peptides on Antigen Presenting 
Cells (APC). Receptors on the membrane of the immune cells are responsible for sensing and translating information 
from the extracellular matrix into the cell. Upon antigen recognition the immune cell orchestrates a spatio-temporal 
organization of its membrane bound proteins that form a protein interface, known as the Immunological Synapse (IS) 
[7] . More broadly, intercellular signaling in a functional IS relates to the formation of large protein domains [H [3] , 
whereas the time scale for their formation and the cluster-to-cluster interaction plays an important role in determining 
the overall cell signaling mechanism. This in turn depends on characterizing the dynamics of the pattern itself, which 
requires us to consider cellular membrane deformations by the receptor-ligand interaction and active cytoskeleton 
processes, leading to fluid motion in the synaptic cleft and thence IS dynamics [8]. 

In the widely studied T-cells, the compartmentalization of membrane-bound protein patterns into different protein 
domains on the cellular scale leads to the formation of Supra Molecular Activation Clusters (SMACs) [HE]. In 
particular, T-Cell Receptors (TCR) form bonds with the peptide Molecular HistoComplex (pMHC) on the APC, while 
Leukocyte-Function-Associated antigen-1 (LFA)-integrin on the T-cell bind with Intercellular Adhesion Molecules 
(ICAM) |3]. A few seconds (0(1 s)) |9] after membrane-to-membrane contact sub micron protein clusters are formed 
that start to translocate (0(1 min)) |3]. This is followed by long range transport and a concomitant coarsening to 
form large-scale protein domains at longer times (0(40 min)) [HilS]- Observations of the T-cell IS show a central 
accumulation of TCR-pMHC, surrounded by a donut-shaped preferential protein domain of LFA-ICAM [3]. 

Understanding the biophysical basis for protein patterning by deciphering the quantitative rules for their formation 
and motion [15] is a first step in characterizing recognition and communication in the immune system. A particularly 
interesting question in this regard is the role of passive physicochemical processes relative to active motor-driven 
processes in generating these patterns m- Recent experiments suggest that early on during the process, active 
processes may not be important, and it is only later that the protein pattern in the T-cell membrane is subject to 
cytoskeletally generated centripetal transport HMj. The question of characterizing the mechanics of the IS patterns 
has led to range of mathematical models that take one of two forms: those that treat the system as a collection of 
discrete units [25H29] or as a continuum |30l [321 [33]. While these models are capable of explaining the spatial 
patterning seen in the IS, they all neglect the fluid flow in the synaptic cleft and thus rely on ad-hoc assumptions for 
the characteristic time scales over which the patterns form, and use approaches based on gradient descent [301 1331433] 
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or stochastic variations of energy minimization of the membrane coupled to protein kinetics |25H28j . 

Here, we provide a description of the passive responses in the IS which includes the mechanical forces due to 
stretching and bending of the cell membrane which are driven by protein attachment and fluid flow, which itself 
causes flow of the trans-membrane proteins. This requires that we integrate cell membrane bending and tension, 
viscous flow in the membrane gap and protein attachment-detachment kinetics, and allows us to capture the essential 
spatiotemporal protein dynamics (nucleation, translation and coalescence of protein clusters) during the formation of 
SMACs. Furthermore, we show that our description of the passive dynamics in the IS implies that the slow dynamics 
of fluid flow can limit the rate of protein patterning, without evoking any active cytoskeletal processes. 


MATHEMATICAL MODEL 
Membrane mechanics 

In Fig. [2 we illustrate the interaction between a T-cell and an antigen seeded bilayer, which mimics the most 
commonly used experimental setup laiiniiiiiiis], and describes the components in the mathematical model described 
below. Once the T-cell is close to the bilayer (Fig. I) the membrane-bound receptors form adhesive bonds with their 
ligand counterparts in the bilayer, which pull the membranes together and squeezing the fluid out of the cleft. 

When these two types of receptors form bonds with ligands, they get compressed or stretched. We assume that their 
spring stiffnesses are inversely proportional to the protein length li that may vary among different protein 

types m- The subscript i = 1 corresponds to the TCR-pMHC complex and i = 2 corresponds to the LFA-ICAM 
complex. Ci = Ci{x, y, t) is the number of attached proteins per surface area (associated with at the total equilibrium 
receptor density Co), their deformation creates a local pressure ^ Ci{x,y,t)Ki{li — h). This pressure deforms the cell 
membrane, approximated here as a bilayer with a bending stiffness Bm = i 2 (f^T 2 )i with E the Young’s modulus, b 
the membrane thickness and v the Poisson ratio (see Supplementary Information (SI)), and a mechanical response 
quantified by 
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p{x,y,t) = + KCi{h - h) + K-^C 2 {h - h) (I) 

h 

where p is the pressure difference across the membrane, and h = h{x^ y, t) is the height of the fluid-filled synaptic cleft. 
By scaling the membrane gap with the longest protein bond l 2 , the lateral lengths with the cell size L and p with 
the spring pressure C 0 KI 2 yields two non-dimensional numbers; B = i2(i-!^2^)kCoL^ describes the relative importance 
of pressure generated by membrane bending and the protein spring pressure, and li/h = 15nm/45nm = 1/3 is the 
ratio of the natural length of the proteins. We focus here on the limit when membrane bending dominates, but we 
show in the SI that the influence of membrane tension smooths some of the small scale pattern features. Active 
cytoskeletal forces would appear as additional source terms in p, but has been neglected below as we focus on the 
passive dynamics. 


Hydrodynamics 


Any membrane deformation initiates fluid motion and give rise to hydrodynamic forces in the synaptic cleft, which 
consequently affects the membrane dynamics. In typical experiments, the synaptic pattern has a lateral size L 
comparable to the cell size (« 10 pm), while the cleft has a height comparable to size of the longest protein bond 
{I 2 = 45 nm). Thus the aspect ratio of the IS is small I 2 IL ^ 1. When combined with the fact that at these small 
length scales, the flow in the synaptic cleft is viscously dominated, we may use lubrication theory |36| to simplify the 
equations governing fluid flow. Under the assumption of a local Poiseuille flow |36j in the membrane gap assuming 
no-slip at both surfaces, where the bilayer of the upper cell membrane can deform but the bilayer on the supported 
glass plate is immobilized. This leads to a single non-linear scalar partial differential equation for the thin film height 
h{x, y, t) [37] similar to that used in other elastohydrodynamic phenomena |351 - HU| 


dh ^ f ^ 

^ = v.(t^lv 

dt \ 12fi 


,i.e. 




h + kCi{Ii — /i) + ^“(^2(^2 — 

h 


( 2 ) 
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(a)Sketch of the interaction between a T-lymphocyte cell (T-cell) and an glass supported antigen 

seeded bilayer. 



^8 nni 


(b) Close view of the synaptic cleft formed between the T-cell membrane and the glass supported bilayer. 


FIG. 1: (a) Sketch of the interaction between a T-lymphocyte cell (T-cell) and a supported antigen-seeded bilayer. The cell 
size L ~ 10 fim and the two membranes are separated by transmembrane receptors bound to ligands in the bilayer, (b) 
Close-up schematic view of the synaptic cleft formed between the T-cell membrane and the glass supported bilayer. The cell 
membrane has a thickness 6 « 8 nm and the membrane gap height is given hy h — h{x,y,t). The trans-membrane receptors 
form bonds with the ligands in the bilayer and have lengths and concentrations, TCR — pMHC « 15 nm, C\(x,y,t), and 
LFA — ICAM ~ 45 nm, C2{x,y,t). During protein bond formation and depletion, the cell membrane deforms generating a 
viscous flow u(a:, y, z) in the membrane gap. The flow generates a viscous frictional force F), parallel with the glass supported 
bilayer that acts onto the cell membrane and the transmembrane proteins and thus affects their motion. Any deformation of 
the membrane generates a restoring elastic bending force Fb, while the deformation of the TCR-pMHC and LFA-ICAM bonds 
generates a spring force Fk- 


where Eq. [^follows by using Eq. [^for pressure (p), where /r is the fluid viscosity. Given the lack of evidence for 
water permeation across the membrane we neglect this effect, as well as thermal fluctuations of the membrane since 
these will be strongly damped out by enthalpic protein binding. 


Protein kinetics 


We only follow the dynamics of the membrane-bound proteins that can bind and unbind from their complementary 
ligands, which is equivalent to stating that the number of these proteins involved in the binding kinetics is large 
compared to the free proteins in the cytoplasm. In the membrane we assume the total number of membrane-bound 
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proteins per unit area to be constant and given by where i = 1 is corresponding to TCR and i = 2 is corresponding 
to LFA. Of these, the number density of bound receptors is denoted by Ci{x, y, t), which can diffuse and get dragged 
due to the fluid flow in the membrane gap, or be actively transported by the cytoskeleton. Their dynamics can be 
described mathematically by a reaction-convection-diffusion equation, which accounts for these effects in addition to 
the binding and detachment of proteins, and in dimensional form reads 

^ = ^VP-VC, + V- (D,Va+’^^^{C,Vh{h-h))]+{C,,o-a)K°'^{k)-aK‘^ffik). (3) 

dt y \ y J 

The first term on the right side is an advective term due to the fluid flow in the synaptic cleft driven by local pressure 
gradients associated with membrane deformation. The flow generates a Stokesian drag on the proteins proportional 
to their size. The second term is a membrane protein flux due to molecular diffusion D^VCi, where the diffusion 
coefficient Di = {li/li)D is assumed to be inversely proportional to the protein length following the Stokes-Einstein 
equation. Alternatively, the membrane diffusivity can be influenced by the membrane anchors, but our results are 
fairly insensitive to the molecular diffusion term (see SI) and we ignore them here. The third term on the right side 
is a drift in response to membrane deformation at a rate {CiS/h{h — U) [301 [Mli where k^T is the thermal 

energy. The last two terms correspond to receptor binding at a rate ~ Ci)K°'^ and unbinding at a rate 
The kinetic rates ATf" and are described in terms of the mean first passage time over an energy barrier mma, 

with a distribution centered around the natural protein length {li) and being a function of li/l 2 — h, given by 

Bond formation: = — exp 

Bond depletion: {li) = = -exp 


where is the kinetic time. To favor protein binding for h ^ li, we assume that proteins lose their bonds three 
times slower (Sr^) [37] than the rate at which they form. Although the exact form of these rates are not known, 
if we assume that the off-rate increases with spring tension, so that proteins would unbind as h li and h ^ U 
and in its simplest form given by a constant off-rate {<Joff = oo) in Eq. (see SI). Experiments show that the the 
different protein pairs form non-overlapping patterns [2ll3ll20j. which we mimic in the choice of the width of the kinetic 
distributions a on = 0.2 and <7off = 0.6 |35j . By narrowing the distributions generate wider protein free areas that 
separate TCR-pMHC and LFA-ICAM rich regions. In contrast, increasing the distribution widths make the different 
protein species overlap, which is unrealistic. We focus here on protein transport due to physicochemical processes 
driven by protein binding and membrane deformation and have neglected the role of active cytoskeleton dynamics in 
the cell e.g. polarized release of T-cell-receptor-enriched microvesicles [24], endocytosis and exocytosis [5|. 




DIMENSIONAL ANALYSIS AND SCALING LAWS 
Dimensional parameters 

The material properties of the cell, the fluid and the proteins that are relevant to the IS and needed as input into 
Eq. 1-4 are summarized in Table [l| as reported in previous work in the literature. 


Dimensionless numbers 

It is natural to scale the horizontal length scales using the cell size, i.e. [cc, y] ^ L, the height of the synaptic cleft 
using the typical protein length i.e. h ^ I2, and the pressure by the local receptor force/area, i.e. p ^ C0KI2 = Po, 
and time by a viscous time, i.e. 

In Eq. [I]|^ the use of the scaled variables p{x,y,t) = p*{x,y,t)pQ = p*{x,y,t)CoKl2, h{x,y,t) = h{x,y,t)*l2, 
X = x*L, y = y*L, t = t*T^, Ci{x,y,t) = (^x,y,t)*Co yields six non-dimensional numbers that govern the dynamics 
of protein patterning, as shown in Table II They are: B = which describes the ratio of pressure generated 

by membrane bending and the protein spring pressure, I1/I2 is the relative ratio between the natural length of the 










5 


Description 

Notation 

Reference 

Fluid viscosity 

/i = 4 X 10““^ Pa ■ s 


Cell membrane Young’s modulus 

E = (0.08 - 80) X 10®Pa 


Membrane thickness 

6 = 8 X 10® m 


Poisson ratio 

V = 0.5 

m 

Bending modulus 

= T2(f^ = 4.5 X (10“®i - 10“!®) 

J [30l|32] 

[S2] 

Protein stiffness (Hookean spring) 

K = 1.2 X 10“®iV/m 

[iHiisn] 

m 

Equilibrium number density TCR 

Cpo = Co = 2 X 10^^m“® 

[3] 

Equilibrium number density LFA 

C 2 ,o = 2 X Co = 4 X 10^'‘m“® 

m 

Natural TCR-pMHC length 

h ~ 15 nm 

m 

Natural LFA-ICAM length 

I 2 = 45 nm 

m 

Membrane protein diffusion coefficient 

D — 5 X 10“^®m/s® 

[501EI] 

Kinetic on-rate 

Tk = ri = T2 — 1.1 X (10“® — 10“’^)s 


Kinetic off-rate 

'Toff = r^ff =rk/3s 

m 

Cell diameter 

L = 10 pm 

m 

Hydrodynamic time scale 

r.= c^=S.7x 10“® s 


Thermal energy 

kbT = 4.34 X 10“®V 


Distribution width on-rate 

^on ~ 0.2 


Distribution width off-rate 

0 

II 

o' 

b 


Pressure scaling 

Po = C 0 K.I 2 = 10.8 Pa 



TABLE I: Description of the material parameters that appear in Eq. 1-4. 


Description 

Membrane bending/protein stretching 
Aspect ratio membrane height/length 
Protein aspect ratio TCR-pMHC/LFA-ICAM 
Diffusive/advective time scale 
Protein sliding mobility/protein diffusion 
Hydrodynamic/kinetic time scale 


Non-dimensional number 

D BfP z 3 

^ ~ kCqL^ ^ 

^ = 4.5 X 10“^ 
lx — lx = 1. 

h h 3 

Pe = = 5 X 10^ 

M = = 2.0 

Dfi 

^=^ = Lrc^=3x(10-3 


2 X (10“’’ 


- 10 ) 


10 “®) 


TABLE II: By substituting the scaled variables in Eq. 1|4 p{x,y,t) = p*{x,y,t)po = p*{x,y,t)CoKl 2 , h{x,y,t) = h{x,y,t)*l 2 , 
X = x*L, y = y*L, t = Ci{x,y,t) = Ci{x,y,t)*Co gives the non-dimensional numbers above. These non-dimensional 

numbers characterize the relative influence of membrane mechanics, protein kinetics, geometry and hydrodynamics. 


proteins, I 2 IL is the aspect ratio of the membrane gap, Pe = ^ is the ratio between advection and diffusion, 

M = is the ratio between protein diffusion and protein sliding mobility, t = ^ — ratio between 

the local hydrodynamic time and the kinetic time (Table |^). As we will show, our results are insensitive to 
variations in Pe, M and initial conditions (see SI), and only the dimensionless numbers B and r control the qualitative 
aspects of our phase space of patterns. 

The two important dimensionless numbers B and r can described the potential variations in the membrane proper¬ 
ties and/or the protein biochemistry across different experiments. In particular, the membrane properties depends on 
its composition, where the presence of inclusions e.g. cholesterol, peptides, proteins, can alter its stiffness, r depends 
on the fluid in the synaptic cleft and the biochemistry of protein binding. In particular, if t > I bonds form rapidly 
relative to the time for fluid flow in the cleft which is then rate limiting, and conversely when r < I, fluid flow is fast 
relative to bond formation which is then rate limiting. 


Length scales 


Two characteristic lengths are observed in the IS, the micro-cluster scale 4 and the large domain scale L. From Eq. 
[^we derive a scaling law for the cluster size, by balancing the spring pressure and bending pressure Bmh/lc^ ~ ^ 0^2 
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t= 5 min U) min 17 min 5 min 10 min 17 min 



AdoOted /rxxn Kaizuka et al., PNAS, 2007 


FIG. 2: Comparison between an experimental (left) and numerical (right) realization of the TCR-pMHC and LFA-ICAM protein 
patterning in the IS. The simulations are based on Eq. |H4| allowing fluid flux at the edge of the IS, where the height and number 
of proteins per membrane area is fixed. The protein pattern is determined by varying B and r. In the experiment a T-cell 
interacts with an antigen seeded lipid bilayer [20]. The upper row shows the density of bonded TCR-pMHC, the middle row the 
bonded LFA-ICAM proteins and the last row their merged RGB-channel. The right panel shows the numerical simulation with 
B = 2 X 10“® and r = 15 where the dimensional times correspond to the non-dimensional times (^)® x t = [1.7, 3.3, 5.7]. All 
other non-dimensional numbers are reported in Table [Hj The white circle illustrates the numerical domain and the experiment 
and simulation are compared at the same instances in time. At short-times, protein clusters nucleate on the membrane, with a 
dynamics given by the interplay between membrane mechanics, protein kinetics and fluid flow. At late times protein clusters 
interact and coalesce into large spatial patterns that mimic pSMAC and cSMAC structures. A ’’donut shaped" LFA ring 
surrounds a dense TCR region at the center of the synapse at late times. 


that leads to 


/ \ I 


(5) 


For the simulated Bm (SI) the deformation length varies between = 70 — 200 nm i.e. in dimensionless 

units I* = si « 0.02 — 0.06 for B G [10“®, 10“^], qualitatively consistent with experimental observations I21j . 


Time scales 


Protein patterning at the micro-cluster (1^) size occurs on short time scales (tc), while patterning at the cell scale 
(L) occurs on long time scales (r^). Fluid continuity and force balance embodied in Eq. [^yields a short time scale 
Tc corresponding to drainage on the micro-cluster scale Ic, given by 


= 12(f) 

h 


= 12 ( 


Bm , 
CokI^' 


CokIo 


Substituting in parameter values yields Tj, « 0.1 — 1 s i.e. in dimensionless time units r* = 12(f )^ 
SI). Fluid drainage on the cellular scale L yields a long time scale given by 


( 6 ) 

24 — 240 (see 


= i2(^; 

h 


= 12 (-) 2 _(^ 


Substituting parameter values yields tl ~ 40 min i.e. in dimensionless units t£ 


5 X 10^. 


( 7 ) 
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t = 2 min 



4 min 



12 min 


: 

>/ 


i<': 




22 min 




^ ^ ^ . 
/ / ^ /'TTA 

/• / J M 


¥|:. 


43 min 


\ \ \ i i t * 

. X X L ^ ^ . 

• X * f r ^ • 

- \ ^*t! . 

■ ^ '*' ' 
f fT\ \ \ \ y 


xl0^2 Pa 

r 




(b) (c) 

FIG. 3: (a)-(c) Simulation of Eq. |l|4| with B = 2x 10“® and r = 15 and the other dimensionless numbers are reported in Table 
2. (a) Contour plots of the time history of the pressure (p) and the velocity {h?Pp), superimposed with the velocity vectors. 
The second row shows the corresponding protein pattern of TCR-pMHC and LFA-ICAM, see Fig. 3 for color scale. At short 
times (t < 4 min) the nucleation and coalescence of protein domains generates a local flow field. At late times (t > 12mm) 
a global centripetal flow is generated that ” compress “ the TCR cluster radially generating a ”bulls-eye”-like protein pattern, 
which becomes unstable at t « 60min. (b-c)The total number of attached receptors (b) TCR-pMHC and (c) LFA-ICAM. (c) 
Direct comparison between of the total number of attached TCR in the IS in simulation and in the experiment by [3] shows 
that the results are in good agreement for t < 20 min. This suggests a ’’waiting-time" for the active cytoskeleton processes 
for protein domain organization and that passive dynamics suffices to describe the short-time formation and organization of 
protein domains. But the long-time IS dynamics and its stability is suggested to be dominated by active processes. 


NUMERICAL EXPERIMENTS 


Methods 


To solve the nonlinear system of Eq. T|4 we numerically discretize these with a finite element method (see SI) 


in two-dimensions, which gives the membrane topography in three-dimensions. For consistency with experimental 
observations, the simulations are performed in a circular domain that capture the central region of the cell-to-cell 
contact, which is assumed not to be influenced by the motion of the cell leading edge. At the edge of the IS the 
membrane is assumed to be torque free with no bending moment (V^/i = 0) and at a constant pressure (p = 0), 
which allows fluid flux through the boundary. The membrane is pinned at the edge {h = O. 5 / 2 ) and the equilibrium 
number of proteins per membrane area at that given height (Ci = C 2 = O.OlCo) see [35] and SI for details. The 
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FIG. 4: Phase space that characterizes the different regimes of membrane protein patterns as a function of B = 

(see Table 2), here presented on logarithmic axes. The simulations are based on Eq. ll4 and the 


and r = — = — T , 

patterning is measured at t = 40 min where a synaptic pattern is typically formed in experiments n El US II y, i.e. in 
dimensionless units (L/hYt = 16. Two different protein patterns are identified; large diffuse patches and dispersed kinetic 
clusters, which are categorized into three regimes. In the diffusional dominated limit (r < 0.3) large diffusive patches are 
predicted that translocate on the membrane. A transition to a dispersed protein pattern is observed for r > 0.3. In the 
intermediate regime (0.3 < t < 3), long-lived LFA clusters form on the membrane. When the protein dynamics is an active 
process (r > 3) micro-scale TCR clusters nucleate and coalesce as they are transported radially forming a central dense pattern. 
In the kinetic regime we see that the cluster size varies as a function of B, similar to our scaling prediction ~ . Note that 

at equilibrium, all simulations predict a flat membrane with a single protein phase for the case where fluid flux at the edge of 
the IS is free and the membrane height and number of proteins per membrane area is fixed. 


membrane is initialized with six small Gaussian shaped bumps of different widths (« O.IL) and amplitude ((0.075 — 
0 . 1 ) 12 )- Additional information about the numerical method |48j . |49j . parameter sensitivity and alternative boundary 
conditions are in the SI. 


RESULTS 

Within the phase space of B and r, we start by considering a cell that has a stiffness that scales with the thermal 
energy « kgT and binding rates that are similar to those reported in experiments m « 10 with an 

association constant « 0.1M“^ giving « 10“®. We note that the hydrodynamic time scale is larger « 3 x 10“^s 
than Tfe suggesting that the IS dynamics is rate limited by the fluid flow i.e. r ^ 1, which we verify below. 

In Fig. we show the time evolution of the IS for these parameters {B = 2 x 10“®, r = 15) and we see that 
the qualitative behavior of our model is consistent with the observed asymmetric IS dynamics seen multiple times 
[3 nn [m na Ha (see Movie 1), and recovers the temporal dynamics and the cluster sizes seen in experiments, 
associated with the presence of dense non-overlapping regions of TCR-pMHC and LFA-ICAM, which vary with time. 
At short times dispersed micron-sized protein clusters nucleate on the membrane, with a characteristic cluster size 
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« 1/zm (containing « 160 proteins). These protein clusters are transported by the centripetal fluid flow generated by 
membrane deformation. At long times, we see the appearance of larger spatial protein domains, with a ” donut-shaped “ 
LFA-ICAM structure (peripheral SMAC) surrounding a dense central domain of TCR-pMHC (central SMAC) (Fig. 
[^. This similarity is particularly striking since we did not evoke any active processes that are present in a cell. 

The two-dimensional simulations of the trans-membrane proteins allow for a direct comparison with the asymmetric 
IS found experimentally H mi [13 Hi sg. To illustrate how these transport processes are correlated with domain 
coarsening, we show the pressure and velocity fields in Fig. At short times (f < 4 min) the nucleation and 
coalescence of protein domains at a length scale « Ic generates a local flow field, while at long times {t > 12 min) the 
flow occurs over a global length scale « L wherein the centripetal flow moves the clusters to the center of the domain 
and coarsens the protein pattern. In Fig. Wp we directly compare the dynamics of the TCR clusters in the simulation 
with experiments |3]. With increasing time, the number of attached TCR rapidly increases upon first contact as micro 
clusters nucleate. A distinct peak in the number of attached TCR is observed around t « 5 min in Fig. [^, followed 
by a decay in the number of attached receptors over longer times. The agreement with experiments for t < 20 min 
is striking since no active processes are evoked and suggests that the slow dynamics of fluid drainage in the synaptic 
cleft limits the rate of protein patterning during the early stages of IS dynamics. 

At longer times (t > 20 min) the results of the simulation and experiments deviate from each other, indicating an 
important role for active processes to stabilize the dynamical synapse. Over this period (« 60 min), a distinctive 
feature in the experiment [3j is the appearance of a stable dense circular region of TCR-pMHC surrounded by a 
’’donut-shaped” ring of LFA-ICAM. Compared to the TCR, the attached LFA display a different dynamics as they 
increase monotonically in time (Fig. [^) and around t « 60 min saturates the nearly flat membrane. A similar time 
evolution is also observed in the experiment by |3| , but their choice of scaling makes a direct comparison challenging. 

Moving beyond the direct comparison with experiments, we turn to a qualitative phase-space of protein patterning 
characterized by T,B,Pe,M, initial conditions and boundary conditions. Our simulations show that the pattern 
dynamics are insensitive to variations in Pe, M and the initial conditions (SI). However, the scaled membrane stiffness 
(B) and the ratio of time seales (t) are the main players responsible for variations in the patterns. In Fig. we 
show this in terms of a phase diagram of pattern possibilities illustrated by snapshots of the protein distributions at 
t = 40 min, a stage corresponding to a mature IS nuns HD]. 

Two distinct protein patterns may be identified corresponding to either large diffuse domains or a dispersed micro 
cluster phase. We can further categorize the latter into two distinct regimes. For r < 0.1 the membrane proteins 
fail to form an IS and their dynamics are primarily dominated by diffusive fluxes and the results are insensitive to 
B. For T > 0.3 islands of non-overlapping micro-scale protein clusters form different shapes on the membrane. For 
0.3 < T < 3 long-lived LFA clusters form at the center and at the edge of the membrane. In this regime, kinetic 
processes and diffusive fluxes make comparable contributions. By further decreasing the kinetic rate (r > 3) the 
protein dynamics become hydrodynamically limited with a sharper protein interface. In this regime, a large central 
domain of TCR with a few internalized LFA micro-clusters form on the membrane, which is surrounded by LFA alike 
the IS. We emphasize that at very long times the equilibrium state corresponds to a nearly flat membrane adhesively 
bound by either TCR or LFA to the bilayer. However, a change in boundary condition that replaces the constant 
pressure along the edge with a vanishing fluid flux, i.e. Vp • n = 0 where n is normal vector at the boundary, leads to 
an arrested inhomogeneous protein pattern (see SI and Movie 2). However, this late stage regime does not influence 
the initial nucleation and growth of protein domains. 

Our calculations of the protein patterns show that the formation of a synapse-like protein pattern only occurs in 
the hydrodynamically limited regime for r > 0.3. In this regime, protein clusters nucleate at short-time t « 1 min 
forming a patchy pattern, with a characteristic cluster size that scales as 4 ~ [^' These micro scale 

protein clusters move centripetally by the self-generated fluid flow since membrane deformation by protein binding 
displaces the interstitial fluid and generates flow, which assists sorting and formation of protein domains. Cluster 
translocation leads to self-interactions and the formation of large protein domains at long times t « 30 min with 
the characteristic ” donut-shaped “ LFA domain that surrounds a central domain dens in TCR (see Fig. 13, similar in 
structure to what is often referred to as a peripheral-SMAC and a central-SMAC in experiments isdinisEffliii]. 


DISCUSSION 

To get at an accurate description of the spatiotemporal dynamics of protein patterning in the IS we have formu¬ 
lated and solved a minimal mathematical model that account for membrane mechanics, protein binding kinetics and 
hydrodynamics, while setting the stage for the quantification of passive and active mechanisms in the IS. Our theory 
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captures the length and time scales of protein patterning seen in experiments by only accounting for the passive 
processes. Our scaling laws for the size of protein clusters, as well as short and long time protein patterning dynamics 
are corroborated in simulations without ad-hoc physical assumptions. In particular we show that slow dynamics of 
fluid drainage in the synaptic cleft can account for the time scales of protein patterning. Direct comparison of our 
computations with experiments by [3] suggests that at early times passive dynamics suffices to describe the formation 
and organization of trans-membrane receptors, and suggests a natural time scale for when active processes come into 
play. Our passive model of the immune-cell synaptic cleft is a simplification, where we have neglected the mecha¬ 
nisms by which receptor binding generates signaling that triggers internal activity e.g. actomyosin polymerization, 
endo-/exo-cytosis, release of TCR through microvesicles, local recruitment of integrins etc. Since all these effects 
can influence the patterning dynamics, to challenge our passive physicochemical theory and to help identify the key 
biophysical process underlying the formation of the IS, we now turn to some experimentally testable predictions. 

First, a characteristic spatial scale for membrane deformation is predicted by Ic = where Bm is bending 

stiffness, Cq protein number density and k protein stiffness. Since Ic is fairly parameter insensitive, modifying cell 
membrane rigidity (wheat germ agglutinin (WGA) [44] ). the protein number density (corralling [46] 1 or protein 
stiffness (linker length [3^) would only produce moderate changes in cluster size. 

Second, two time scales are derived for the short and long time dynamics. At short time protein clusters nucleate 
Tc « and at long time and length scales large protein domains form tl « CoKh ' 

contrast to the prediction for Ic, both Tc and tl are sensitive to changes in protein number density (Co), protein (k) 
and membrane stiffness {Bm), which can be experimentally changed by corralling, linker-length and WGA and will 
change these three parameters, respectively. Thus, our theory predicts that the time scales for the IS can be changed, 
without much variation in the spatial features. 

Third, our numerical simulations predict only protein domains for < Tk, identifying protein kinetics as a critical 
component in the IS formation. Thus, proteins need to bind faster (r^) than the characteristic fluid flow (r^) to 
form a protein pattern. By changing the adhesion molecules the kinetic rate can be varied and can be modified 
by playing with the protein number density (corralling) or protein stiffness (linker length), whereas fluid viscosity is 
expected to be challenging to alter. 

Fourth, the effective boundary condition at the periphery of the synaptic cleft is found to be a key component 
in the longevity of the pattern. Simulations allowing fluid flux through the edge of the IS show that the SMAGs 
become unstable at long times. The formation of a tyrosine phosphatase network at the synapse periphery generates 
additional resistance to fluid drainage and may limit the rate of mass flux. Thus, the proteins at the boundary of the 
IS are predicted at one component that regulate its stability and suggests that a disruption of this protein network 
would affect its longevity. 

Fifth, the fluid motion in the membrane gap has hitherto not been quantified. Such experiments may be feasible 
with quantum dot tracing techniques m) and may shed new light on the fluid pathway during the patterning. Fluid 
can either become trapped in the inter-membrane gap, internalized by the cell or escape at its edge. Another time scale 
appears for a flow though a porous cell protein network (^ which depends on its permeability Kp [m“^] 

and thickness Im [m]. 

Sixth, we predict nucleation, translation and sorting of protein clusters in the absence of active processes. Recent 
observations by m of non-immune cells show protein patterning and makes an experimental platform ideal to 
challenge our spatiotemporal predictions. 

Our mathematical model presented here is a minimal and general theoretical skeleton for an accurate description of 
this class of cell-to-cell interaction phenomena, which is likely useful beyond the IS and understand broader aspects 
of cell adhesion, communication and motility. 
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SUPPLEMENTARY INFORMATION FOR ELASTOHYDRODYNAMICS AND KINETICS OE PROTEIN 
PATTERNING IN THE IMMUNOLOGICAL SYNAPSE” A. CARLSON AND L. MAHADEVAN 


Description of movies 

Movie 1: The dynamics of protein patterning for the case when the fluid flux at the boundary is free to vary, but 
the pressure is fixed. This causes the pattern to eventually decay. 

Movie 2: The dynamics of protein patterning for the case when the fluid flux at the boundary vanishes. This causes 
the pattern to eventually get arrested. 


Problem parameters 

Table |T] summarizes material properties that are relevant to the IS synapse, as reported in previous work in the 
literature and are used as inputs to Eq. 1-4. 


Boundary conditions 

The boundary condition at the edge of the IS critically affects the final protein pattern since it reflects the role 
of different biophysical processes associated with membrane deformation, fluid flow and the number of proteins per 
membrane area. Three different types of boundary conditions for the membrane edge can be prescribed, as given 
below: 

Pinned: "S/^h = 0, h = constant, Ci = C 2 = constant 
Clamped: Nh ■ n = 0,h = constant, Ci = C 2 = constant 
Free: = V^h • n = 0, VJ, • n = 0, 

where n is the boundary normal. In a slight abuse of notation, we denote Pinned and clamped as corresponding to 
fixing the density of proteins at a given height at the membrane edge, and further letting the torque vanish or fixing 


TABLE I: Description of the material parameters that appear in Eq.(l-4). 


Description 

Notation 

Reference 

Fluid viscosity 

/i = 4 X 10“^ Pa ■ s 


Cell membrane Young’s modulus 

E = (0.08 - 80) X 10®Pa 


Membrane thickness 

6 = 8 X 10® m 


Poisson ratio 

u = 0.5 


Bending modulus 

^ 12(1-1^2) — 4.5 X (10 

10"^®) J [301132] 

m 

Protein stiffness (Hookean spring) 

K = 1.2 X 10~® Y/m 

[291130] 

m 

Equilibrium number density TCR 

Ci,o = Co = 2 X lO^^m-® 

m 

Equilibrium number density LEA 

Ca.o = 2 X Co = 4 X 

m 

Natural TCR-pMHC length 

h ~ 15 nm 

m 

Natural LFA-ICAM length 

I 2 = 45 nm 

m 

Membrane protein diffusion coefficient 

D — 5 X 10“^®m/s® 

[SolEI] 

Kinetic on-rate 

r*, = n = ra = 1.1 x (10“® - 10' 


Kinetic off-rate 

U"// = Toff =Tk/Ss 

m] 

Cell diameter 

1 / = 10 fim 

m 

Hydrodynamic time scale 

T. ~ c:.i 2 - 3-7 X 10-® s 


Thermal energy 

kbT = 4.34 X 10"®V 


Distribution width on-rate 

^on ~ 0.2 


Distribution width off-rate 

aoff = 0.6 


Pressure scaling 

Po = C 0 KI 2 = 10.8 Pa 
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the angle at the edge. For a shear and moment free membrane edge, the edge is free and further we assume that there 
is no protein flux at the boundary. 

Any of these three set of boundary conditions for the membrane height can be prescribed with the additional 
boundary conditions for the fluid motion that read 

Free fluid flux: p = 0 

No fluid flow fu = Oj: Vp • n = 0. 

In the case where there are few proteins at the membrane edge, the pressure is prescribed {free fluid flux) to allow 
fluid flow into or out of the membrane gap. If the IS is sealed off by a dense protein network, it is hard for the fluid 
to escape at the edge and no fluid flux is a natural prescription. 

In the main text, we have used a pinned membrane with a constant pressure at the edge, allowing mass fluid flux 
through the boundary of the domain. The use of these boundary conditions is based on experimental observations, 
where only at late times a protein network surrounds the IS. One implication of using a pinned and free fluid flux 
boundary conditions is that the protein pattern does not stabilize, which can however be arrested by prescribing 
no fluid flow at the free boundary. We note that it is possible to extend the mathematical model to account for a 
free boundary problem for the location of the edge itself, but we avoid this scenario here as it is not relevant to the 
dynamics of the IS and generates a significant numerical complication, as a dynamic mesh is needed to track the 
membrane edge. 


Computational methodology 

The governing equations (Eq. 1-4) were solved with the open-source finite element toolbox femLego [48]. A first 
order semi-implicit Euler scheme is used for time marching and all variables are discretized in space using piecewise 
linear functions. The non-linearity together with the sixth order derivatives in Eq. 2 makes it challenging to solve. 
Therefore, we decompose this into three equations, for the Laplacian of the height (V^h), the pressure (p) and 
the height {h). These three equations are coupled with two additional equations for the proteins, and are solved 
simultaneously using a Newton iteration method [49] . 

We performed a convergence study of the ID results in both time and space; different spatial Ax/L = 
[0.007,0.0025,0.001] and temporal = [100,400,1000] resolution show no noticeable change in the results. The 
2D results presented here use a mesh size {Ax/L = 0.006) and a time step ((^)^At* = 0.02). 


Dependence of dynamics on ratio of hydrodynamic to kinetic time scale r 


To investigate the nature of the spatiotemporal evolution of the trans-membrane proteins we vary r while keeping 
all other parameters fixed (Fig. S.ll. In Fig. S.l we see that when r ^ 1 the patterns are kinetically limited. 
Large protein patches nucleate on the membrane that slowly drift by diffusion. In contrast, increasing the role of 
hydrodynamics by the increase in r leads to a patchy protein pattern of receptor micro-clusters that are separated 
by a sharp interface. The clusters move centripetally (see Fig. S.l), causing the pattern to coarsen as they coalesce, 
which leads to the formation of large protein domains. At equilibrium the membrane is nearly flat and saturated by 
a single protein species. 


Parameter sensitivity - diffusion, sliding, initial condition, off-rates and membrane tension 


For given initial conditions, protein diffusion, sliding and advection can influence the dynamics of patterning. To 
quantify the influence of these properties on the resulting protein patterns, we separately turn off these effects. In 
Fig. S.2 we show the results when protein sliding is turned off {M~^ = 0), in Fig. 8.2] we show the results when 
protein advection is turned off (4hi/L^ = 0), and in Fig. 8.2: we show the results when protein diffusion is turned off 
{Pe~^ = 0) (Fig. 8.21). What is clear from Fig. 8.2 is that none of these parameters has any significant contributions 


in the kinetic regime (r = 3.0), where macroscopic patterns persist. 

In order to determine the influence of diffusion {D) on the resulting dynamics, we varied D over two orders of 
magnitude (Fig. l8.2^ -g). By changing the diffusion coefficient the value of both Pe and M change. Increasing D 
makes the interface between the two boundaries more diffuse, while decreasing D makes the boundaries sharper. 
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FIG. S.l: Time history of the attached proteins and the membrane topography as a function of r = [0.03,0.3,3.0,30]. r = 
^ ^ relative importance between the local viscous time (r^) and the kinetic time rate (ta,). B = 2 x 10“® 

and the other non-dimensional numbers are given in Table 1. The simulations are based on Eq. 1-4. The color-scale for the 
density of bonded LFA (green) and TCR (red) is shown to the upper left corner and the scale bar for the membrane height 
(black-white) is shown in the upper right corner. For r ^ 1 the dynamics are hydrodynamically limited and no protein clusters 
are predicted. In contrast, for r > 0.3 clusters of TCR and LFA nucleate at short-time and translocate centripetally at long 
times forming large protein domains. 


Comparing the results for D = 5.0 x 10“^^m^/s and D = 5.0 x I0“^^m^/s we notice that besides this quantitative 
feature, the results are indistinguishable and diffusion does not strongly influence the protein patterns over this range. 

To quantify the influence of the initial conditions, we perform simulations with three different initial conditions 
(Fig. |S.2| r-j). Initially the membrane has six small Gaussian shaped bumps of different widths (« O.IL), with an 
amplitude ((0.075 — O.I)Z 2 )- In the sub-figure to the lower left in Fig. 7 the bumps on the membrane are inverted 
compared to the simulation to the lower right. The result presented in the middle sub-figure shows a simulation 
result with an initial membrane shape with six Gaussian bumps at different positions than shown in the left and 
right sub-figure. Although the detailed shape of the pattern is slightly influenced by the initial condition, the overall 
dynamics is robust to these changes. 

We have in this work assumed that the kinetic binding and unbinding rates are described by the means passage 
time over an energy barrier Eq. 4, which leads to an Gaussian distribution for the on/off rates centered around li. 
The off-rates may be a function of the tension in the proteins with a probability of unbinding that increases with 
the tension up to a given threshold. Eq. 4 generates an effective kinetic rate that takes the form of a double-well, 
while an off-rate based on the tension in the proteins would remove the two minima and the probability of unbinding 
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B = 2x 10-'^ 

= 0.51 
r = 3.0 

Pe-^ = 1.85 X 10“® 


All contributions Sliding Advection Diffusion 



M ^ = 0 _ Pe ^ = 0 


L2 


Influence of diffusion coefficient: 

B = 2x 10-'^ 

T = 3.0 




Pe-^ = 1.85 X 10-^ Pe-^ = 1.85 x 10“® Pe”^ = 1.85 x 10“® 
M-i = 5.1 M-i=0.51 M-i = 0.051 


Influence of initial condition: 

B = 2x 10“® 
r = 3.0 

Pe-^ = 1.85 X 10"^ 



FIG. S.2: Influence of protein diffusion, sliding and advection and the initial condition on the predicted numerical results at 
time t = 23 min i.e. in dimensionless units = 7.0. 


approaches a constant value as the proteins are further stretched/compressed. The simplest form of a tension based 
off-rate is to let Uo// = oo, where the effective rate (Ff°" — K°^^) becomes a shift of the gaussian for the on-rate 
and the probability of unbinding becomes constant for large protein deformation. We have performed additional 
simulations to verify that our results are not very sensitive to the from of the off-rate, which is demonstrated in the 
second row in Fig. S.3. Although the detailed shape of the pattern is slightly different, the overall dynamics is robust 
predicted in the simulation. 

Since the membrane has a fluid-like nature, there can also be an influence in the pressure from membrane tension 
and an additional term enters into Eq. 1 

21 

p{x,y,t) = + KCi(h — h) + K-^C 2 {h — I 2 ) (S. 8 ) 

where 7 is the membrane tension [N/m\. In the tension dominated limit the length scale for membrane deformation 
1 

scales as • Scaling pressure with the characteristic spring pressure KI 2 C 0 yields another dimensionless number 

r = ^ which is the ratio between pressure from membrane stretching and the pressure from deforming the 

protein springs. In Fig. S.3 row 3-5 we demonstrate the influence of membrane tension by varying F = [10® — 10®] 
e.g. 7 = [4 — 0.04] X 10“®7V/m in dimensional units. For B = 2 x 10“^ we note that as F < 10“® the spatiotemporal 
dynamics is dominated by membrane bending (Fig. S.3). If the membrane tension is increased, larger protein domains 
appear « Fs and if the membrane becomes too stiff the pressure generated by the protein springs is not sufficient to 
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deform the membrane F > 10 ^ 





FIG. S.3: Time history of the attached proteins and the membrane topography as a function of off-rate (Joff = oo (row 2) and 
membrane tension F = [10^ — 10®] (row 3-5) for i? = 2 x 10“^, r = 3.0, Pe = 5 x 10^ and M = 2.0. F = is the ratio of 
pressure from the membrane tension and the protein spring pressure. The color-scale for the density of bonded LFA (green) 
and TCR (red) is shown to the upper left corner and the scale bar for the membrane height (black-white) is shown in the upper 
right corner. These snapshots in time correspond to the dimensionless times (^)^ x t* = [14, 28, 71,142]. 
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FIG. S.4: Influence of boundary conditions on the IS dynamics (r = 3.0, B = 2 x 10 ®) at times (a) t = 10 min = 2.8), 

(b) t = 23 min = 7.0), (c) t = 47 min ((^)^t* = 17)and (d) t = 230 min ((^)^t* = 69). At the edge the membrane 

moves freely using a shear and moment free boundary condition, with no fluid flow and a no-flux boundary condition is 
prescribed for the TCR-pMHC and LFA-ICAM proteins. Contrary to the pinned membrane (Fig. S.ll, which allow in- and 
out-fluid flow, the protein pattern is arrested at long times. 


Sensitivity to boundary conditions 

To illustrate how the boundary condition can affects the simulation results in Fig. |S.4| we show a sequence of 
snapshots of a simulation with a membrane that is allowed to move freely at the edge {shear and moment free) with 
no-flux of proteins and no fluid flow. Comparing these results with the case when proteins are free to diffuse through 
the boundary, it is clear that the boundary condition affects the protein patterning and serves to arrest the protein 
pattern at long times. 



